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ABSTRACT 

We determine the most likely values of the free parameters of an N-body model for the 
Galaxy developed by Fux via a discrete-discrete comparison with the positions on the 
sky and line-of-sight velocities of an unbiased, homogeneous sample of OH/IR stars. 
Via Monte-Carlo simulation, we find the plausibilities of the best-fitting models, as 
well as the errors on the determined values. The parameters that are constrained best 
by these projected data are the total mass of the model and the viewing angle of the 
central Bar, although the distribution of the latter has multiple maxima. The other 
two free parameters, the size of the Bar and the (azimuthal) velocity of the Sun, are 
less well-constrained. The best model has a viewing angle of ~ 44°, semi-major axis 
of 2.5 kpc (corotation radius 4.5 kpc, pattern speed 46 kms -1 kpc -1 ), a bar mass of 
1.7x 10 10 M and a tangential velocity of the local standard of rest of 171 kms -1 . We 
argue that the lower values that are commonly found from stellar data for the viewing 
angle (~25°) arise when too few coordinates are available, when the longitude range is 
too narrow or when low latitudes are excluded from the fit. The new constraints on the 
viewing angle of the galactic Bar from stellar line-of-sight velocities decrease further 
the ability of the Bar's distribution to account for the observed micro-lensing optical 
depth toward Baade's window : our model reproduces only half the observed value. 
The signal of triaxiality diminishes quickly with increasing latitude, fading within 
approximately one scaleheight (SS3°). This suggests that Baade's window is not a very 
appropriate region to sample Bar properties. 
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1 INTRODUCTION 

Evidence has accumulated over the last five years that there 
is a triaxial structure in the inner Galaxy (see Gerhard 1996 
for a review) . The study of the galactic Bar received a large 
stimulus when the COBE-DIRBE data and the derived 
models became available (Dwek et al. 1995). Earlier, star 
counts, gas dynamics and three-dimensional stellar kine- 
matics had been analysed. Although most studies agree on 
the presence of a Bar and roughly on its orientation, the 
exact viewing angle, size and shape remain a matter of de- 
bate. In this article we use for the first time a global set of 
stellar line-of-sight velocities in the galactic Plane to deter- 
mine the values of these parameters, by comparing the set 
to a triaxial N-body model of the Galaxy (Fux 1997). 

N-body models are important for the study of the 
dynamics of the triaxial Galaxy, because they are self- 
consistent and have known formation- and evolution his- 



tory. This is opposite to self-consistent Schwarzschild-type 
models (Schwarzschild 1979), for which we know only the 
present and the future. Schwarzschild-type models can be 
aimed directly at fitting observations. True N-body models 
can much less easily be 'steered' that way and comparing 
them with observations is difficult. This holds especially for 
the Galaxy, where stellar kinematic data are always discrete 
and one is therefore faced with the problem of comparing 
two distributions of discrete data points. One could smooth 
and normalize the (projected) N-body model if it has suf- 
ficient particles and sample this probability distribution at 
the observed points. The result is the probability of the data 
given the model. One can also smooth the observations to 
obtain a velocity profile, as is often done with observations 
in Baade's window, and compare that to the corresponding, 
smoothed, N-body profile. Statistical tests then yield the 
probability that the two distributions are the same. 

We have the opportunity to use a new stellar-kinematic 
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data set, homogeneous and unbiased, with highly accurate 
on-the-sky positions and line-of-sight velocities (Seven- 
ster et al. 1997a,b; S97A, S97B). This data set (hereafter 
AOSP (Australia telescope Ohir Survey of the Plane)) is 
particularly suited for constraining dynamical models for the 
Galaxy, not only because of its high accuracy but also be- 
cause of the intrinsic properties of the stars. The AOSP sam- 
ple consists of so-called OH/IR stars, observable throughout 
the Galaxy. These are evolved, intermediate-mass stars and 
their distribution, spatial as well as kinematic, traces closely 
the global galactic potential (Habing 1993, Frogel 1988) and 
is relatively relaxed. OH/IR stars have circum-stellar en- 
velopes due to mass loss and the outflow velocity of those 
circum-stellar envelopes can be obtained directly from the 
spectra. The outflow velocity is related to the stellar mass 
and thus to the age of the star, in a statistical sense (van 
der Veen 1989; see Sevenster 1997). This allows, for example, 
determination of the changes in the dynamical distribution 
with time. 

The average surface density of the sample is of the order 
of one star per square degree. This means that two neigh- 
bouring stars cannot be assumed to sample the same velocity 
profile, which is implicitly required to smooth the data. Also, 
it is not necessary to smooth the model completely, because 
we want to determine the probability of the model given the 
data, rather than the other way round, or the probability 
that model and data have the same distribution. We use 
a method to scale an N-body model (Fux 1997) to match 
the data, via an implementation of a direct discrete-discrete 
comparison (Saha 1998). The model was chosen from a range 
of N-body models because it reproduces best the combina- 
tion of the COBE-DIRBE surface-density map (in the K 
band (2.2//m)) and other observations (eg. the local disper- 
sions and density, Fux 1997). It is therefore most representa- 
tive for the AOSP sample, because this comes from exactly 
the same intermediate-mass, evolved stellar population that 
dominates the near infra-red surface density observed by 
COBE. 

In Section 2 we describe briefly the general method for 
the Galaxy model-data comparison, in Section 3 we describe 
the detailed implementation for the given data and model we 
use. The results are presented and discussed in Sections 4,5. 
In Section 6 we calculate the implications for gravitational 
micro lensing toward the galactic Bulge and we finish with 
conclusions in Section 7. 



2 THE METHOD 

2.1 Determining the best fit 

To compare the six-dimensional N-body model (cartesian 
coordinates x,y, z,u,v,w) and the three-dimensional data 
(galactic longitude I and latitude b, and line-of-sight veloc- 
ity V), the model is projected according to : 

x = x cos 4> + y sin (f> y — y cos <f> — x sin <f> + Rq 



u cos (f> + v sin < 



£ = arctan 
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} — Vq sin £ cos b. 



(1) 



The four scaling parameters <j), Rq, /v ,Vq are the free pa- 
rameters of the model. The viewing angle (f> is the orien- 
tation of the Bar with respect to the line of sight to the 
galactic Centre (if = 0° the Bar points toward the Sun). 
Rq enters equation (1) as if it were the distance of the Sun 
to the galactic Centre, but it determines the size of the 
Bar (£, b become smaller with larger _Rq ) . If in the initial 
model the semi-major axis is a, then in the scaled model it 
is (a x 8kpc/_R ). This is the consequence of having only 
(£, b, V) to fit; if we had more coordinates there would be an 
extra free parameter /r (the true size scale) and Rq would 
be the true distance to the galactic Centre. In equation (1) 
this parameter /r is hidden, as it enters the numerator as 
well as the denominator for all three coordinates and thus 
is of no consequence. The velocity factor /v determines the 
total mass, M p oc / v (8kpc/i?0) (maintaining virial equilib- 
rium). Vq is the azimuthal velocity of the Local Standard 
of Rest (LSR) . The four free parameters are all weakly cor- 
related, as is clear from the fact that we have only three 
quantities to fit. In equation (1) we see immediately that 
the spatial distribution of the data places no constraints on 
/v and Vq. The line-of-sight velocities do constrain R Q 
which enters the coordinate-transformation term as well as 
the solar-motion-correction term in the expression for V. 
The viewing angle <f> plays a prominent role in all terms of 
equation (1). 

We determine the values of the free parameters that op- 
timize the model-data fit with a so-called genetic program- 
ming method (Charbonneau 1995). We divide the three- 
dimensional data-space in B = NeNbNy boxes and deter- 
mine the number of model- and data particles in each box 
i, rrii and Oi, respectively. The joint probability W that the 
data (in total O particles) and the model (M particles), 
projected on the data-space, arise from the same underly- 
ing distribution function, is given by the following formula 
(Saha 1998) : 



(rm + Oi)'. 
(rriiW.) 



C 



M\Q\(B - 1)! 
(M + + B-l)\ 



Oi,rrii = Oi,m;(0, i? Q , V©, /v) 



(2) 



Note that this equation is symmetric in the model- and 
data-terms. W is robust against outliers in the data (or in 
the model if M < O which is unlikely ever to be the case) 
and also against unphysical solutions, such as putting all M 
model particles in the box with highest Oi (Saha 1998). 

If vtii = or Oi — 0, no contribution to the likelihood 
W is made; the term within the product in equation (2) 
equals one. Preferably, B < M and B < O so that we have 
as few boxes as possible without information content. On 
the other hand, we want to prevent any smoothing of the 
data, so that Oi&l for all i and B > O. In other words, 
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B ~ O, which results in B < M as M » O. A funda- 
mental constraint on the box size comes from the demand 
that, within a box, the distribution function, that M and O 
derive from, is constant. This constraint is much harder to 
quantify in practice, because we do not know the distribu- 
tion function. In general it is also in conflict with the first 
constraint. If we make the boxes so small that the distribu- 
tion function is truly constant within each of them, they will 
not all contain at least one star. A proper balance has to be 
found between the two constraints. From tests we find that 
results are robust over a large range of B - 20 < B < 80. 
For both larger and smaller B/O the plausibilities start go- 
ing to zero, although the best-fit values remain constant for 
smaller B/O. Judging by the values, biases and plausibilities 
of our tests, B/O ~ 4 is optimal. 

In principle all model particles could be used for opti- 
mizing equation (2) , but in practice we use random subsets 
with M — 16384, saving calculation time without losing pre- 
cision. We carry out the optimization for 100 such subsets 
to obtain a number distribution of best-fit values for each 
parameter (see Fig. 1) . We then determine the medians of 
these distributions, or all local maxima if there are more 
than one (only <f> as will be clear later), thus finding one or 
more best fits (<j>, Rq, Vq, fv). 

2.2 Determining the plausibility 

Having found a best fit to the data, we want to know whether 
it is also a good fit, within the limits of the model. Via 
Monte-Carlo simulation (eg. Press et al. 1986) we determine 
the "intrinsic" W distribution for the best-fit model. Sets 
of O model particles are randomly selected from the (entire) 
best-fit model and the corresponding value of W is calcu- 
lated. The resulting distribution thus gives the W values 
for the case we know that the model and the "data" are 
the same. We then determine the percentage P(W) of these 
values that is lower than the W of the real data. A high 
value of P(W) means that the result is significant; the fit is 
as good as can be expected for that particular model. We 
carry out this loop (the left branch in Fig. 1) 25 times, with 
50 different random subsets each time, to obtain (P). We 
will call (P) the plausibility of the fit. Roughly, models with 
(P) < 10% are not acceptable, models with (P) > 50% are 
optimal. 

2.3 Determining the errors 

To quantify the errors that are connected with the determi- 
nation of a best-fit model, we again draw random samples 
of O particles from the (entire) best-fit model and find the 
best fit for these fake data sets. The mean value for each of 
the free parameters from these fits, and the dispersion in the 
values, show the intrinsic accuracy of the fitting procedure. 
For each quantity A, we thus find the la error cta and differ- 
ence between the best-fit value and the Monte-Carlo mean, 
in terms of the la error, the bias Xa = (A mean — A& t )/aA- 
In practice, we use 25 independent samples, created in such 
a way that O n M = . 

3 IMPLEMENTATION 

3.1 Model and data particles 

The N-body model we use is the model "m08t3200" (at 3.2 
Gyr in the simulation) developed by Fux (1997). It contains 
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Figure 1 Schematic view of the procedure to determine the 
W— distribution (left branch) and the intrinsic error on the free 
parameters A (right branch). The double arrows indicate the 
steps that involve optimizing (of W). First the maximum of 
W(ip, Rq, Vq, fv) (equation (2)) is determined with a genetic 
algorithm for 100 different subsets of the model. From the re- 
sults, best-fit values A are determined (see Section 2.1). From 
the model, scaled with parameter values A, 50 random samples 
of O particles are drawn and the value of W for each of those 
"data sets" given the scaled model is calculated. The whole left 
branch is executed 25 times, yielding the plausibility (P) (see 
Section 2.2). On the other branch, 25 random samples of O par- 
ticles are drawn from the scaled model. For each of those, again 
W(4>, Rq, Vq, fv) is optimized, ie. using the initial model. The 
resulting distributions in </>,Rq, Vq,/v yield the errors <r in and 
the biases X on the fitted values. 

a bar that formed spontaneously (without imposed triaxial 
potential) from an axisymmetric distribution of stars, in a 
disk and a spheroid, plus a dark halo. The stellar part con- 
sists of 100,000 particles of 6.57 x 10 5 M . Corotation is at 
5.4 kpc (determined from the moments of inertia) and the 
bar's semi-major axis is 3 kpc (out to the start of spiral 
arms, around the inner ultra-harmonic resonance) in the 
initial model. The circular velocity, at the radial range in 
which it is constant, is 218 kms -1 , so the pattern speed of 
the bar is of the order of 40 km s _1 kpc -1 . The in-plane axis 
ratio of the model Bar is approximately 0.5. Symmetry with 
respect to the plane was imposed during the simulation. The 
ratio of the corotation radius to the semi-major axis is 1.8, 
which together with an exponential density profile, makes 
this a late-type Bar (Elmegreen 1996, Noguchi 1996). 

The AOSP data set consists of 507 OH/IR stars (Sec- 
tion 1) with measured on-the-sky positions, accurate to 0'.'5, 
and line-of-sight velocities, accurate to 1 kms -1 . The ob- 
servational errors are effectively zero in this analysis and will 
be neglected in the rest of this article. The properties of the 
stars in the AOSP sample (Section 1) allow us to gain more 
information from this modelling than just the best fit to the 
full data set. The outflow velocity of the circum-stellar enve- 
lope, of the order of 15 km s _1 , is roughly proportional to lu- 
minosity (van der Veen 1989), mass and age (since these are 
all the same parameter, given a certain stellar-evolutional 
phase). The relation should be applied in a statistical sense. 
Stars with higher outflow velocities can hence be detected 
out to larger distances (^12 kpc), on average, than those 
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with lower outflow velocities (^10 kpc, Sevenster 1997), pro- 
vided they are observed with the same flux-density cut-off 
(S97A,S97B). The spectra of some of the stars (<20%) show 
only one of the usual two peaks, thus not allowing for a de- 
termination of the outflow velocity and an accurate line-of- 
sight velocity (S97A, S97B). For these stars, as a group, the 
velocity accuracy is of the order of the average outflow veloc- 
ity, 15 kms -1 . We applied the method to the total AOSP 
sample as well as to subsamples, to see how their different 
properties influence the fit (Section 3.3). The line-of-sight 
velocities of the stars are given with respect to the local 
standard of rest (LSR, see S97A). Throughout this article, 
we will use V© and "solar motion" to indicate the azimuthal 
motion of the LSR. 

3.2 The final runs 

The outflow velocity of OH/IR stars correlates, as mentioned 
above, with average distance. We used that property by run- 
ning the fitting program with various subsets of the data. 
In Table 1 the windows on the sky, the number of boxes in 
each coordinate, the number of stars O and the total num- 
ber of boxes B that we used in the various runs are given. 
For all runs, M = 16384 so that the fits do not depend on 
the number of model particles within the window. 

As mentioned in Section 2.1, optimally B/O ~ 4. The 
"standard" run bd uses the total AOSP sample and the 
total region on the sky. Setting (N(, AT b , iVv)=(20,6,15) was 
found optimal in tests for this run. They correspond to bin 
sizes of 2 ° 75, 1 ° 1 and 40 km s _1 . For the other runs, we try 
to stick to those values of B/O and the bin sizes as closely 
as possible given the individual run's characteristics. 

In run b, we use the largest possible longitude range 
symmetric about the galactic Centre. In run bdd, we use 
only sources with double-peaked spectra, because they have 
the best-defined line-of-sight velocities, with, for this appli- 
cation, negligible errors (1 km/s, S97A, S97B). 

In run bdl we use only sources with outflow velocities 
between lkms -1 and 15kms _1 and in run bdh only those 
with outflow velocities larger than 13kms _1 . The ranges in 
outflow velocities overlap slightly in order to retain sufficient 
stars for the comparison. 

The bdh contains stars at large distances (Section 3.1) 
and the fit for this run will be sensitive to the full morphol- 
ogy of the Bar. For run bdl, the data particles and model 
particles probably do not trace the same distances; the sub- 
sample is complete out to ^10 kpc (Section 3.1). We there- 
fore have run bdll with a distance window (0 kpc - 10 kpc) 
for the model. 

To assess the degree to which the fit is determined by 
the stellar velocities or by their positions, we ran the pro- 
gram with increased velocity- (vel) and spatial resolution 
(mor), respectively. Those two runs have the same B and 
O to facilitate comparison. We also "switched off" the kine- 
matics completely (morx), as an extreme test, keeping in 
mind that Ny — 1 does not satisfy the demands on the box 
sizes. Finally, to judge better the suitability of the AOSP 
sample for constraining the free parameters of the N-body 
model, in ptl we used a different sample of OH/IR stars (te 
Lintel Hekkert et al. 1991). This sample is incomplete in a 
number of ways. First, it is an IRAS-selected sample, which 
means it is incomplete at very low latitudes where confu- 



sion limited the number of point sources the IRAS satellite 
detected. Second, the velocity coverage changes with longi- 
tude. Third, the infra-red selection (see te Lintel Hekkert 
et al. 1991) introduces a inhomogeneous distance sampling 
that is difficult to quantify. The reason for using this sam- 
ple nevertheless, or actually because of all this, will become 
clear. 

In fixf we fix all parameters, except <f>, at the values 
found by Fux (1997) for the model m08t3200; the value used 
for V© is that of the local circular velocity in his best-fit 
model. Furthermore, in frxl we fix <f> at 25°, in frx.2 at 45° 
and in in fix3 at 65°, all three maxima in the distributions 
of best-fit values of 4> in various runs. 

3.3 The search ranges 

In preliminary tests, we found that the optimal search ranges 
are (0° - 90°) for <f>, (6 kpc - 10 kpc) for R , (160 kms^ 1 - 
230 kms" 1 ) for V Q and (0.25 - 2.25) for M p . Since scaling 
the total mass means scaling the potential and hence the 
velocities squared, in practice we determine this parameter 
by scaling a velocity factor /v between (0.5 - 1.5), M p oc /y. 
The only exception is run morx, where we use (0.3 - 3.0) for 
/v- For the observable quantities, these ranges span amply 
the likely real values. 

4 RESULTS 

Figures 2-5 show the results of the 100 W-optimizations 
(see Fig. 1) for some of the runs. The values of the free pa- 
rameters are found to be virtually uncorrelated, so that we 
can determine the maxima or medians of these distributions 
for each parameter separately to obtain the best-fit mod- 
els (4>, Rq, Vq, /v). Clearly, the distributions of R & and V© 
do not have well-defined maxima (inside the search range). 
We use the median value for those parameters, as well as 
for /y. The latter's distribution has the best-defined max- 
imum of all four parameters (except in morx) but owing 
to its symmetric and smooth distribution the maximum is 
the same as the median value in all cases. We find several 
clear local maxima only for (j> and determine these from the 
derivatives of the cumulative number distributions (ie. the 
unbinned actual number distributions; see the dashed curves 
in Figs 2-5). 

In Table 2, we give the sets of best-fit values for each of 
the runs described, along with the plausibility (P) of this 
best fit and the spread in the plausibility ap . For each of the 
parameters, we also give the la error, <ja, and the difference 
between the best-fit value and the bias Xa, determined as 
described in Section 2. It should be noted that for best-fit 
values at an upper or lower boundary of a search range, Xa 
is necessarily large, negative or positive, respectively. In run 
morx without the kinematics, the parameters V© and /v 
(in italics) are not constrained at all. 

In general, it is clear that Rq and V© are not optimally 
constrained by our method and/or data. Their best-fit dis- 
tributions do not show clearly isolated maxima and from 
Table 2 we see that these parameters have average biases of 
~0.5, as opposed to ~0.2 for <f> and /v. In terms of the ra- 
tio of (7 a to the corresponding search range, fy is very well 
confined with at on average 5% of the search range. For the 
three other parameters this ratio is 20-25%. There is good 
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Table 1. Windows for the runs (M = 16384) 
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Figure 2 The cumulative number distributions from the 100 
best-fit values for the four free parameters from run bdd. For 
<j> (panel a), the dashed curve gives the (arbitrarily normalized) 
derivative of the cumulative number distribution and hence shows 
the (unbinned) actual number distribution. 
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Figure 3 As Fig. 2, for run bdl. 



agreement to within ~ lcr between the various models (ex- 
cept b, bdll, morx and ptl, this will be discussed later) on 
the values of the parameters Rq (8.9-9.6), Vq (164-179) and 
fv (0.90-0.95). For <j> the situation is considerably different. 
Interestingly, despite the fact that the best-fit distribution of 
4> mostly has several maxima, the subsequent Monte-Carlo 
analysis of each of the solutions shows that some are remark- 
ably well confined. In Fig. 6 we show the histograms of the 
values for <j> in Table 2 as well as of all the values for <j) occur- 
ring in the 100 optimizations for the six bd* runs. The three 
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Figure 5 As Fig. 2, for run ptl. 

peaks around 25°, 45° and 65°obviously instigated the runs 
fixl— 3. We disregard the fourth peak at 85° as this really 
indicates an axisymmetric solution. For <^>=90° (side-on) the 
surface-density profile is completely symmetric. 

In principle run bd, using the most datapoints and the 
largest longitude window, should give the best results. Of 
its three solutions, the 44° fit coincides most closely with 
fits from several other runs (Fig. 6). Most noteably, the runs 
with increased velocity accuracy, either from the data (bdd) 
or from the method (vel), have a solution for similar <j> with 
very small bias in <j>- Also bdh, that samples completely all 
distances throughout the Bar, has a fit with 0=45°. These 
four <j> ~44° runs give similar values for R and V©, and 
reasonable agreement for fv- Moreover, the bdh 0=45° fit 
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Table 2. Best-fit parameter values (local maxima for <j>, medians for Rq, M p , Vq) and plausibilities 
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has very high plausibility, as has fix2 with respect to fixl, 3. 
The 0=44° result for bd (to be called bdl) therefore gives 
the best scaling parameters for Fux's N-body model. 

4.1 The degeneracy in <f> 

Obviously, there are other values of cj> that give equally rea- 
sonable fits (judging by the errors, biases and plausibilities) . 
In the lower panel of Fig. 6 we see that 4> even follows an in- 
trinsically quadrumodal distribution ! It was shown by Zhao 
(1997b) that, from surface density only, the viewing angle 
cannot be determined uniquely. Different density models for 
the Bar can give exactly the same projected densities with 
different viewing angles. It is expected, however, that no 
such degeneracy would exist when optimizing a given spa- 
tial density model, or when including global kinematics con- 
straining self-consistent model dynamics. The fact that we 
still find a degeneracy is probably at least partly dictated by 
the limited extent of our longitude window at positive longi- 
tudes. 'Asymmetry profiles' (positive-to-negative-longitudc 
ratios of surface density or mean velocity etc.) are crucial in 
constraining cj> and these profiles are cut short. In the next 
paragraph we discuss what influence this may have. 



In runs ptl and bdll R© is ~2cr below the value of the 
other runs. The samples used in both runs have larger (ap- 
parent) scaleheights than the whole population of evolved 
stars; in ptl because the sample is incomplete in the plane 
and in bdll because the distribution of older stars simply has 
larger (intrinsic) scaleheight. Run bdl, however, the equiva- 
lent of bdll but without integration limit (see Section 3.1), 
results in an average value for R Q . Possibly, this can be 
explained by the facts that we are not using the integra- 
tion limit and that the observed longitude range extends to 
+10°, ie. the near tip of the Bar is outside the window, as 
we mentioned before. This has a two-fold effect. First, the 
surface-density asymmetry profile of the cut-off bar may 
look like the profile of an entire bar with larger viewing an- 
gle. We suggest that this is the origin of the <j> values around 
25°, as these are seen in bdl, morx and ptl. All samples 
that contain the bdl stars are clearly incomplete at larger 
distances to some extent, most of all ptl (the integration 
limit can not be estimated for this sample) . Apparently, the 
kinematics suppress the tendency for <j) to be ~25° (see also 
Fig. 6) if the incompleteness is not too strong, so that we 
only see it in run morx, that does not include kinematics, 
and in bdl and ptl with severe incompleteness (not taken 
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Figure 6 In the top panel, the histogram of all best-fit values 
of <j> (from Table 2) is shown; in the bottom panel of all 100 
W— optimizing values from the six bd* runs. Note the narrow 
peak at 45° , the wide peak around 65° and the small peak around 
25°, as well as the absence of values below 20°. The peak at 85° 
probably indicates near-symmetry (90°) of some data sets. 

into account. Second, one can imagine that with the far end 
of the Bar 'cut off' in the data by the flux density limit and 
the near end by the longitude limit, the result may look like 
a smaller bar, in other words : larger R Q . This may cause 
the different R© in bdl and bdll. One would expect the 
larger R Q in ptl if it didn't have the longitude extend up 
to +25°. 

Just as bdll and ptl have larger than average scale- 
height, bdh has smaller scaleheight. Indeed, in bdh the 
value for R is relatively large. The deviation is not as large 
as for bdl, as the bdh sample spans the whole Bar and puts 
stringent constraints on the surface-density asymmetry pro- 
file. 

The values for Vq are rather low compared to the 200±20 
kms -1 currently accepted (see Dehnen & Binney 1998; 
Feast & Whitelock 1997; Rohlfs et al. 1986). However, these 
determinations all assume V is the local circular speed 
which is unlikely, amongst other reasons because the disk 
may be slightly elliptical (Kuijken & Tremaine 1994). These 
authors advocate an average local circular velocity of 200 
km s _1 and a tangential velocity for the LSR of 180 km s _1 , 
making the low fitted values for Vq more acceptable. The 
difference with the model's circular velocity of 207 kms -1 
is still large. In run b V© is particularly ill-determined, be- 
cause the sini term in the correction for the line-of-sight 



velocities (equation (1)) covers a small range only. 

Very high velocities (up to 450kms _1 ) are present in 
the initial N-body model that were never found for OH/IR 
stars. Altogether, less than ten stars are known at absolute 
velocities higher than 300 kms~ (Baud et al. 1975; van 
Langevelde et al. 1992; S97A). These high model velocities 
are not used in the comparison, however, because they are 
outside the windows defined in Table 1. Those windows are 
the true limits of the observations (S97A, S97B) so increas- 
ing the velocity window would be meaningless. The data- 
model comparison, in particular the determination of fy, is 
therefore not based on the extreme tails, but on the wings of 
the distribution of the bulk of the velocities. The resulting 
/v<l shows that the total mass M p of the initial N-body 
model is somewhat too large. 

4.2 Best model 

As we have argued, the best values for the free parameters 
are given by the 44° solution of bd : bdl. From the free 
parameters we can derive some more interesting properties 
of the best-fitting model. Corotation radius is at 4.5 kpc and 
the semi-major axis of the Bar is 2.5 kpc (cf. Section 3.1). 
The pattern speed is 46 kms -1 kpc -1 for the rescaled local 
circular velocity of 207 kms -1 . Finally, the mass of the Bar 
is 1.7xl0 10 M . The plausibility of bdl is 15±7% , but 
when determining (P) with only the double-peaked stars 
(bdd) we get 68±8% . The single-peaked stars in bd, that 
possibly do not (all) belong to the OH/IR (asymptotic-giant 
branch) star population, apparently decrease the goodness 
of the fit. 



5 DISCUSSION 

The value of 44° we obtain for <f> is large compared to some 
other estimates (16°, Binney et al. 1991 (gas dynamics); 
20°, Binney, Gerhard & Spergel 1997 (integrated light); 24°, 
Nikolaev & Weinberg 1997 (star counts); 20°-30°, Stanek et 
al. 1997 (flux differences between positive and negative lon- 
gitudes)). It should be noted that the lowest value for the 
viewing angle, 16° (Binney et al. 1991), is based on the pos- 
sibly wrong assumption that the CO "parallelogram" (Bally 
ct al. 1988) is formed by gas on the inner cusped orbit (see 
Sevenster 1997). The parallelogram may well be the result of 
a tilt in the inner gas disk (Liszt & Burton 1978). Our high 
(f> is compatible with the COBE E2&E3 models (<f> ~40°) by 
Dwek et al. (1995) and with the value of 35° determined by 
Weiner & Sellwood (1996; from the gas dynamics). Also, Un- 
avane & Gilmore (1998) find from near-infrared starcounts 
that models with viewing angles between 20°-45° are ac- 
ceptable. A much older determination of <j> from the gas 
dynamics even yields 45° (Peters 1975). To explain the two 
local star streams, Kalnajs (1997) argues that the viewing 
angle should be 45° as well. 

Not only the high value for viewing angle, but also the 
value for the corotation radius, or pattern speed, can be 
reconciled with the observed gas kinematics. Various at- 
tempts to model HI and/or CO kinematics in the inner 4 
kpc have given values for the pattern speed ranging from 
fip^igkms^kpc- 1 (Wada et al. 1996) to 63 km s" 1 kpc" 1 
(Binney et al. 1991) or even 118kms _1 kpc -1 (Yuan 1984). 
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Figure 7 Surface density, mean galactocentric line— of— sight ve- 
locity and dispersion (all at 6 = 0°) for model bdl in open cir- 
cles and for data (total OH/IR sample) in filled squares. The 
points shown are obtained via 3-D-adaptive-kernel smoothing. 
Note that hence the distributions are not reliable at the borders 
of the observed region (££7°; ££-42°). 

Hydrodynamic models for the HI longitude-velocity dia- 
gram of the whole Galaxy are illustrated in Mulder & Liem 
(1986), who themselves give as best model parameter 0=20° 
and Rcr ~R© (their Fig. 5), but their model with 0=40° 
and Rqr ~0.5 Rq (their Fig. 9) gives a similarly good and 
arguably better fit to the 3-kpc arm and the central CO 
kinematics (Bally et al. 1988). 

In Fig. 7 we illustrate how the data compare to the bdl 
model. The global agreement is good and also the region 
around I = —20° which may be in the corotation region of 
the Galaxy (Sevenster 1997) shows very similar features in 
all three quantities in model and data. 



Fux himself (1997) derives a best-fit value for <f> of 25° 
with the same N-body model (his values for the other pa- 
rameters are as in fixf, Table 2). We find that the values 
for <f> found in run fixf are no different from the other runs. 
To mimic the model optimization from only the COBE K- 
band surface density (Fux 1997; his value fy comes from 
scaling the velocities to fit the line-of-sight dispersion to- 
wards Baade's window) we introduced run morx where the 
kinematics are 'switched off' completely by setting Ny=l. 
Indeed, one of the solutions gives 0=29° (and R©=9.0 kpc), 
but it carries zero plausibility. As argued before, the global 
stellar kinematics seem essential in determining the viewing 
angle as they provide the necessary constraint to suppress 
the degenerating influence of limited windows and distance 
coverage. In Fux's (1997) case the longitude window (+30°, 
—30°) is large enough to prevent the problems we discussed 
in Section 4.1. However, latitudes |6| <3° are excluded from 
the optimization because the COBE data cannot be cor- 
rected reliably for extinction in the plane. This reminds us 
of run ptl with its incompleteness in the plane and subse- 
quent low viewing angles and incapability of constraining 
the model ((P}=0%). We would argue that many of the 
4> ~25°-30° results found in the literature suffer from simi- 
lar problems. 

The measure for the relative residual defined by Fux, 
R 2 Np . x=300 , is 1.5% for bdl (cf. 0.47% for his best fit for 
m08t3200). 

It has proven virtually impossible to distinguish be- 
tween triaxial and axisymmetric distributions studying the 
distribution of the line-of-sight velocities only (eg. Ibata 
& Gilmore 1995; Dejonghe et al. 1997). We applied the 
Kolmogorov-Smirnov test as used by Ibata & Gilmore 

(1995) , as well as the intrinsically more powerful distance- 
velocity statistic described by Dejonghe et al. (1997), to the 
N-body model. Both only give significant results for very 
low latitudes (|6| < 2°) and even then only with samples 
of at least 1000 stars. The results are very dependent on 
the viewing angle. It is therefore no surprise that Ibata & 
Gilmore (1995) found no significant evidence for triaxiality 
from their intermediate-latitude velocity profiles and that 
neither statistic gives a signal when applied to the AOSP 
sample. However, even though the relatively low surface den- 
sity of the AOSP sample inhibits the construction of velocity 
profiles, it contains sufficient kinematic evidence of triaxial- 
ity, as indicated by the plausibility of >50% . The low lati- 
tudes, very homogeneous sampling and simultaneous fitting 
of the spatial and the kinematical distribution are essential. 

In an earlier stage of this project, we applied the same 
procedure to the Schwarzschild-type N-body model by Zhao 

(1996) in its initial state before evolution. No significant fit 
was obtained for this unmixed model. Fux's N-body model 
has proven to have a solid physical basis and shows many ob- 
served features in a variety of data (Fux 1997; Fux & Friedli 
1996). Also the N-body bar's formation, spontaneously via 
instability of the underlying disk, is one of the probable 
ways to form bars (eg. Sellwood & Wilkinson 1993). The 
fact that for the ptl sample, selected especially for its in- 
completeness (see Section 3.2), no significant fits can be 
obtained, gives extra credibility to the model (as well as 
the method). This all provides proof that m08t3200, Fux's 

(1997) best model, gives a very good representation of the 
six-dimensional Galaxy. 
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Clearly, there is room for improvement, especially to 
lift the ^-degeneracy once and for all. The most impor- 
tant features of a stellar data set used to achieve this are 
a large longitude range and homogeneous sampling of low 
latitudes. More dimensions to limit the degrees of freedom 
are preferred over more objects. 



6 MICRO-LENSING OPTICAL DEPTH IN 
THE LINE OF SIGHT TO THE BAR 

In this section, we want to discuss briefly the micro-lensing 
properties of the scaled model. The observed micro-lensing 
optical depth is as yet unaccounted for by any density model 
for the central Galaxy so it is important to assess this optical 
depth for our best model. In a forthcoming paper, we will 
do this in more detail and also calculate the event-duration 
distribution, that may give insight in the nature of the miss- 
ing optical depth, even though the stellar mass function is 
not known in the model. 

The micro-lensing optical depth r is the probability 
to detect a micro-lensing event at a given instant. From a 
density model, one can calculate the distribution of r on 
the sky; the micro-lensing-optical-depth map. Comparison 
with the measured values (r = 3. 9t\\t x 10 -6 for red clump 
giants toward (£ = 2? 55, b = -3° 64); r = 2.1±° B 4 x 1(T 6 
for main-sequence stars, corrected for blending, toward (t 
= 2? 7, b = -4?1), Alcock ct al. 1997; r = Z.Z±\% x 1(T 6 
toward Baade's window {I ~1°, b ~ -4°), Udalski et al. 1994) 
gives important information about the model. For a wide 
range of models that derive from the COBE maps (Dwek et 
al. 1995), one finds that r mo d to be 2<r lower than r b s (Zhao 
& Mao 1996). The missing optical depth thus has to be 
accounted for, within the limits put by other observations, 
by a component not present in those models; either dark or 
sub-stellar particles or an extra density component (eg. a 
thick disk). 

For the calculation of r one needs to take the bright- 
ness of the lensed sources into account. This is some function 
of their distance D s , so that the optical depth in a certain 
direction also depends on D s . Kiraga & Paczynski (1994) 
hence defined Tp oc D 2+213 , where (3 defines the exact de- 
pendence on D a and from the subscript of r one can im- 
mediately see which dependence was used in the theoretical 
calculation. For (3 = 0, the sources are visible out to infinity; 
we get ro <x D 2 and the proportionality is a straightforward 
volume correction. 

For bdl, we calculated the ro map (Fig. 8), as described 
in Fux (1997). It is well known that only a strong and mas- 
sive bar, with a viewing angle smaller than 20° (Zhao & Mao 
1996; Kiraga & Paczynski 1994; Zhao 1997a; Fux 1997 his 
Fig. 13), can fully account for the measured values. Models 
derived from observed surface-density measurements have 
bars too weak, and mostly at too high viewing angle, to be 
able to account for r (eg. Nikolaev & Weinberg 1997). The 
values we obtain in the direction of the measurements are 
ro ~ 1.5 x 1CF 6 , best compared to the 3.9 x 10 -6 value 
for the clump giants, and r_i ~ 0.9 x 10~ 6 , best compared 
to 2.1 x 10~ 6 . Indeed both ro and t_i (both without the 
contribution of the dark component) are too low to explain 
observations by > 2cr, thus confirming the preliminary re- 
sults presented by Valls-Gabaud et al. (1997). 
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Figure 8 Micro-lensing-optical— depth (ro) map, symmetrized 
in latitude, for bdl, without the contribution of the dark par- 
ticles in the simulation (see Fux 1997). The thick contour is for 
ro = 1.0 X 10~ 6 ; the other contours are spaced by a factor of 1.5, 
decreasing with increasing latitude. 

The distribution of the micro-lensing optical depth is 
not very bulge-like; it is dominated by lenses in the disk. 
There is no significant asymmetry between positive and neg- 
ative longitudes, as would be expected when either lenses 
or sources (or both) have a barred distribution (see Evans 
1994). However, the asymmetries in the optical-depth dis- 
tribution from disk lenses and bar lenses, respectively, have 
opposite signs (Evans 1994). Simple tests, with (3 = 0, — 1, 
show that disk lenses skew the distribution toward positive 
longitudes; bar lenses toward negative longitudes. The net 
effect clearly depends on (3, (f> and density parameters. 

In fact, the asymmetric signal is expected to be largest 
for a viewing angle 45° (eg. Evans 1994). In principle, since 
in Fig. 8 we show the ro map (/? = 0), we expect the distri- 
bution to be skewed toward negative longitudes. Additional 
spiral arms, protruding from the ends of the Bar (see Fux 
1997, m08t3200 model), could counteract this asymmetry. 
All in all, the symmetric appearance of Fig. 8 is understand- 
able. 



7 CONCLUSIONS 

We presented values for the free parameters of a self- 
consistent model of the Galaxy, optimized to fit positions 
and velocities of various sets of evolved stars. These stars 
are representative for the global galactic distribution. The 
method is found to be sensitive to incompletenesses and to 
large-scale kinematics. 

The over-all distribution of the stars is fitted well by a 
bar of the global form of that of the N-body model (Fux 
1997) with a semi-major axis of 2.5 kpc, corotation ra- 
dius of 4.5 kpc, an axis ratio of 0.5 and a viewing angle 
of 44° . The value for the viewing angle is high but not in- 
compatible with previous determinations from stellar data 
as well as gas kinematics. The mass contained within this 
Bar is ~ 1.7 x 10 10 Mq, marginally lower than various other 
derivations (Zhao et al. 1996; Blum 1995; Kent 1992). The 
derived pattern speed for the Bar, Q p = V c (Rcr)/Rcr is 46 
kms -1 kpc -1 . For the solar azimuthal velocity a low value of 
171 kms -1 is found. This is much lower than the local cir- 
cular velocity (207 kms -1 ) in the model potential, implying 
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that the Sun (more precise, the local standard of rest) would 
not be on a circular orbit in this particular model. This is in 
agreement with the findings of Kuijken & Tremaine (1994) 
who find a local circular velocity of 200 kms -1 and V© = 
180 kms" 1 . 

We argue that using low-latitude, unbiased, global stel- 
lar kinematics is crucial to determine the viewing angle <f>. 
The commonly found and accepted low values of (f> ~ 25° 
(see Section 5) should be viewed with caution. Our method, 
applied to a variety of data sets with 'known flaws' or to the 
stellar positions only, shows that indeed these favour view- 
ing angles around the lower value. The resulting fits do not 
have high plausibility. 

One of the reasons is that the signal of the Bar dimin- 
ishes quickly with increasing latitude and thus <j> becomes 
ill constrained. The sample of ptl, that is underrepresent- 
ing the plane below ~3° - approximately one scaleheight 
- gives a null-result. This means that either there has not 
been significant bar-induced thickening in the inner Galaxy, 
or the thickening conspires with the distribution becoming 
rounder. We believe that parameters of the galactic Bar can 
not be reliably constrained without data that trace its inner 
one scaleheight. 

This also means that Baade's window may not be an 
appropriate region to sample the Bar's properties (see also 
Sevenster 1997). Regardless of this consideration, with the 
best-model viewing angle of 45° the Bar does not give signif- 
icantly higher values for r than do a variety of axisymmetric 
distributions (Zhao 1997a; Kuijken 1997) and r mod ~ 0.5T obs 
for our best model. We conclude that (provided that the 
value of T bs is beyond suspicion) the origin of the discrep- 
ancy between current bar models and the observed micro- 
lensing optical depth should be sought in a foreground com- 
ponent - eq. a spiral arm or a thick disk - with larger appar- 
ent scaleheight than the Bar. The most convincing argument 
in favour of a small viewing angle for the Bar (Zhao 1997a) 
is thus taken away. 
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